%% plot_CTBCS_phase3D.m
% 绘制 CTBCS 三维耦合系统的相空间（对应论文 Fig. 3a）

clc; clear; close all

%—— 参数设置 —————————————————————————————————
r1        = -10;      % CTBCS 控制参数 r1
b1        = -3;       % CTBCS 移位常数 b1
X0        = [0.3; 0.5; 0.6];  % 初始向量 [x0; y0; z0]
N_total   = 6000;     % 总迭代步数
transient = 1000;     % 丢弃前 transient 步
%————————————————————————————————————————————————

% 预分配有效点
M  = N_total - transient;
XS = zeros(M,3);

% 初始化：前一、当前状态
X_prev = X0;
X_curr = CTBCS3D_step(X_prev, r1, b1);

idx = 0;
for i = 1:N_total
    % 三维 CTBCS 迭代
    X_next = CTBCS3D_step(X_curr, r1, b1);
    % 更新
    X_prev = X_curr;
    X_curr = X_next;
    % 丢弃暂态后存储
    if i > transient
        idx = idx + 1;
        XS(idx, :) = X_curr';
    end
end

% 绘制三维相空间
figure
plot3(XS(:,1), XS(:,2), XS(:,3), '.', 'MarkerSize', 3)
xlabel('x'); ylabel('y'); zlabel('z')
title('(a) CTBCS Phase','FontWeight','normal')
set(gca,'FontSize',13)
grid on
view(135,25)


%% ———— 子函数 ————

function X2 = CTBCS3D_step(X1, r1, b1)
% CTBCS 三维耦合迭代：F₃ = 3D 量子Logistic，G₃ = 3D‑SIMM
    F3 = quantum3D(X1, 3.99, 30);
    G3 = SIMM3D_vec(X1, 1, 2*pi, 11.5);
    X2 = cos(pi*(r1*F3 + (1-r1)*G3) - b1);
end

function Y = quantum3D(X, r, b)
% 3D 量子 Logistic（论文 Eq.(1) 实数化简化版）
    x = X(1); y = X(2); z = X(3);
    Y = zeros(3,1);
    Y(1) =  r*(x - x^2)       - r*y;
    Y(2) = -y*exp(-2*b) + exp(-b)*r*( (2 - 2*x)*y - 2*x*z );
    Y(3) = -z*exp(-2*b) + exp(-b)*r*( 2*(1 - x*z)*z - 2*x*y - x );
end

function Y = SIMM3D_vec(X, a, b2, c)
% 3D‑SIMM
    x = X(1); y = X(2); z = X(3);
    x1 = sin(b2*z)   * sin(c/x);
    y1 = sin(b2*x1) * sin(c/y);
    z1 = sin(b2*y1) * sin(c/z);
    Y  = a * [x1; y1; z1];
end
